Theo Bakker, lector Learning Technology & Analytics, De HHs
Publicatiedatum
25 juni 2024
1.1 Methode, data en analyse
1.1.1 Toelichting op de methode
Voor de ontwikkeling van prognosemodellen gebruiken we de aanpak van Tidymodels. Tidymodels is een framework voor het bouwen van een prognosemodel. Hiermee verzekeren we ons van een systematische, herhaalbare en schaalbare aanpak.
1.1.2 Toelichting op de data
De basis voor deze analyse is studiedata van De Haagse Hogeschool (De HHs), verrijkt door het lectoraat LTA. De data bevat informatie over de inschrijvingen van studenten in het eerste jaar van de opleiding:
Demografische kenmerken: geslacht, leeftijd, reistijd en SES totaalscore.
Vooropleidingskenmerken: toelaatgevende vooropleiding, studiekeuzeprofiel, gemiddeld eindcijfer in de vooropleiding en eventuele deelname aan het Navitas programma.
Aanmeldingskenmerken: aansluiting (direct na diploma, tussenjaar, switch), dag van aanmelding, aantal parallelle studies aan De HHs en collegejaar.
1.1.3 Toelichting op de analyse
We toetsen in deze analyse Retentie na 1 jaar, voortaan Retentie genoemd.
Retentie is gedefinieerd als ingeschreven staan in dezelfde opleiding in een aansluitend collegejaar. Een wisseling van opleidingsvorm binnen de opleiding, bijv. van voltijd in jaar 1 naar duaal in jaar 2, geldt ook als retentie.
Uitval is het tegenovergestelde van retentie: niet ingeschreven staan in dezelfde opleiding in een aansluitend collegejaar. Een wisseling van opleidingsvorm binnen de opleiding, bijv. van voltijd in jaar 1 naar duaal in jaar 2, geldt niet als uitval.
1.2 Voorbereidingen
1.2.1 Laad de data
We laden een subset in van historische data specifiek voor:
Opleiding: BRV | B Bestuurskunde Overheidsmanagement (BO), duaal, eerstejaars - Retentie na 1 jaar
Toon code
## Laad de data voor de opleidingdfOpleiding_inschrijvingen_base <-get_lta_studyprogram_enrollments_pin(board ="HHs/Inschrijvingen",faculty = faculteit,studyprogram = opleidingsnaam_huidig,studytrack = opleiding,studyform =toupper(opleidingsvorm),range ="eerstejaars")## Herschik de levelsSet_Levels(dfOpleiding_inschrijvingen_base)dfOpleiding_inschrijvingen_base <- dfOpleiding_inschrijvingen_base |>## Maak een eenvoudige succes variabele aanMutate_Retentie(sSucces_model) |>## Maak van de succes variabele een factormutate(SUC_Retentie =as.factor(SUC_Retentie)) |>## Verbijzonder eventueel op basis van het propedeusediploma# Filter_Propedeusediploma(sPropedeusediploma) |>## Maak van de Dubbele studie variabele een Ja/Nee variabelemutate(INS_Dubbele_studie =ifelse(INS_Aantal_inschrijvingen >1, "Ja", "Nee")) |>## Verwijder INS_Aantal_inschrijvingenselect(-INS_Aantal_inschrijvingen) |>## Pas voor een aantal variabelen de levels aanMutate_Levels(c("VOP_Studiekeuzeprofiel_LTA_afkorting","INS_Aansluiting_LTA","VOP_Toelaatgevende_vooropleiding_soort" ),list(lLevels_skp, lLevels_vop, lLevels_vop) )## B Huidtherapie: Filter op uitsluitend studenten met een rangnummer (selectie)if(opleiding =="HDT") { dfOpleiding_inschrijvingen_base <- dfOpleiding_inschrijvingen_base |>filter(!is.na(RNK_Rangnummer)) }
1.2.2 Selecteer en inspecteer de data
We selecteren eerst de relevante variabelen. We verwijderen daarbij variabelen die maar 1 waarde hebben. We bekijken de variabelen in een samenvatting in relatie tot retentie. Daarnaast bekijken we de kwaliteit van de data op missende waarden.
Toon code
lSelect <-c("INS_Student_UUID_opleiding_vorm","CBS_APCG_tf","DEM_Geslacht","DEM_Leeftijd_1_oktober","GIS_Tijd_fiets_OV","INS_Collegejaar","INS_Dagen_tussen_aanmelding_en_1_september","INS_Dubbele_studie","INS_Aansluiting_LTA","INS_Navitas_tf","SES_Deelscore_arbeid","SES_Deelscore_welvaart","SES_Totaalscore","SUC_Retentie","VOP_Gemiddeld_cijfer_cijferlijst","VOP_Gemiddeld_eindcijfer_VO_van_de_hoogste_vooropleiding_voor_het_HO","VOP_Cijfer_CE1_nederlands","VOP_Cijfer_CE1_engels","VOP_Cijfer_CE_proxy_wiskunde","VOP_Cijfer_CE1_natuurkunde","VOP_Studiekeuzeprofiel_LTA_afkorting","VOP_Toelaatgevende_vooropleiding_soort" )## B Huidtherapie: voeg de variabele RNK_Rangnummer toeif(opleiding =="HDT") { lSelect <-c(lSelect, "RNK_Rangnummer")}## Maak een subsetdfOpleiding_inschrijvingen <- dfOpleiding_inschrijvingen_base |>## Selecteer de relevante variabelenselect_at(lSelect) |>## Hernoem variabelen voor beter leesbare namenrename(ID = INS_Student_UUID_opleiding_vorm,Geslacht = DEM_Geslacht,Leeftijd = DEM_Leeftijd_1_oktober,Reistijd = GIS_Tijd_fiets_OV,Dubbele_studie = INS_Dubbele_studie,Collegejaar = INS_Collegejaar,Aanmelding = INS_Dagen_tussen_aanmelding_en_1_september,Aansluiting = INS_Aansluiting_LTA,Navitas = INS_Navitas_tf,APCG = CBS_APCG_tf,SES_Arbeid = SES_Deelscore_arbeid,SES_Welvaart = SES_Deelscore_welvaart,SES_Totaal = SES_Totaalscore, Retentie = SUC_Retentie,Cijfer_SE_VO = VOP_Gemiddeld_cijfer_cijferlijst,Cijfer_CE_VO = VOP_Gemiddeld_eindcijfer_VO_van_de_hoogste_vooropleiding_voor_het_HO,Cijfer_CE_Nederlands = VOP_Cijfer_CE1_nederlands,Cijfer_CE_Engels = VOP_Cijfer_CE1_engels,Cijfer_CE_Wiskunde = VOP_Cijfer_CE_proxy_wiskunde,Cijfer_CE_Natuurkunde = VOP_Cijfer_CE1_natuurkunde,Studiekeuzeprofiel = VOP_Studiekeuzeprofiel_LTA_afkorting,Vooropleiding = VOP_Toelaatgevende_vooropleiding_soort ) |>## Pas CBS_APCG_tf aan naar factormutate(APCG =case_when(APCG ==TRUE~"Ja", APCG ==FALSE~"Nee",.default ="Onbekend")) |>## Geef aan waar missende cijfers in het VO zijnMutate_Cijfers_VO() |>## Verwijder variabelen, waarbij er maar 1 waarde isselect(where(~n_distinct(.) >1)) |>arrange(Collegejaar, ID)## B Huidtherapie: hernoem de variabele RNK_Rangnummerif(opleiding =="HDT") { dfOpleiding_inschrijvingen <- dfOpleiding_inschrijvingen |>rename(Rangnummer = RNK_Rangnummer)} dfOpleiding_inschrijvingen <- dfOpleiding_inschrijvingen |> ltabase::sort_distinct()## Verwijder de basis datasetrm(dfOpleiding_inschrijvingen_base)
Studentkenmerken versus Retentie
Variabele
Retentie
p-value2
Totaal, N = 11361
Ja, N=628 (55%)1
Nee, N=508 (45%)1
Aanmelding
144,37 (63,88)
137,01 (67,09)
0,044*
141,08 (65,41)
Aansluiting
0,008**
2e Studie
0 (0%)
5 (100%)
5 (100%)
Direct
336 (53%)
300 (47%)
636 (100%)
Na CD
21 (78%)
6 (22%)
27 (100%)
Switch extern
190 (57%)
145 (43%)
335 (100%)
Switch intern
29 (57%)
22 (43%)
51 (100%)
Tussenjaar
52 (63%)
30 (37%)
82 (100%)
APCG
0,10
Ja
188 (51%)
178 (49%)
366 (100%)
Nee
427 (58%)
315 (42%)
742 (100%)
Onbekend
13 (46%)
15 (54%)
28 (100%)
Cijfer_CE_Engels
7,11 (1,28)
7,43 (1,06)
0,008**
7,23 (1,21)
Cijfer_CE_Engels_missing
<0,001***
Ja
309 (50%)
307 (50%)
616 (100%)
Nee
319 (61%)
201 (39%)
520 (100%)
Cijfer_CE_Natuurkunde
6,30 (0,90)
5,93 (1,21)
0,48
6,20 (0,98)
Cijfer_CE_Natuurkunde_missing
0,032*
Ja
605 (55%)
500 (45%)
1.105 (100%)
Nee
23 (74%)
8 (26%)
31 (100%)
Cijfer_CE_Nederlands
6,34 (0,92)
6,37 (0,86)
0,57
6,35 (0,90)
Cijfer_CE_Nederlands_missing
<0,001***
Ja
309 (50%)
308 (50%)
617 (100%)
Nee
319 (61%)
200 (39%)
519 (100%)
Cijfer_CE_VO
6,63 (0,46)
6,60 (0,49)
0,16
6,62 (0,47)
Cijfer_CE_VO_missing
0,011*
Ja
246 (51%)
237 (49%)
483 (100%)
Nee
382 (58%)
271 (42%)
653 (100%)
Cijfer_CE_Wiskunde
6,38 (1,14)
6,28 (1,15)
0,30
6,34 (1,14)
Cijfer_CE_Wiskunde_missing
0,002**
Ja
343 (51%)
324 (49%)
667 (100%)
Nee
285 (61%)
184 (39%)
469 (100%)
Cijfer_SE_VO
6,62 (0,46)
6,58 (0,48)
0,12
6,61 (0,47)
Cijfer_SE_VO_missing
<0,001***
Ja
277 (50%)
276 (50%)
553 (100%)
Nee
351 (60%)
232 (40%)
583 (100%)
Dubbele_studie
<0,001***
Ja
1 (7,7%)
12 (92%)
13 (100%)
Nee
627 (56%)
496 (44%)
1.123 (100%)
Geslacht
0,051
M
332 (53%)
298 (47%)
630 (100%)
V
296 (58%)
210 (42%)
506 (100%)
Leeftijd
19,74 (2,58)
19,88 (2,62)
0,32
19,80 (2,60)
Reistijd
44,61 (24,59)
43,36 (25,30)
0,22
44,05 (24,90)
SES_Arbeid
0,00 (0,09)
-0,01 (0,09)
0,11
-0,01 (0,09)
SES_Totaal
-0,03 (0,31)
-0,06 (0,33)
0,092
-0,04 (0,32)
SES_Welvaart
-0,01 (0,14)
-0,03 (0,15)
0,066
-0,02 (0,15)
Studiekeuzeprofiel
0,001**
EM
134 (50%)
133 (50%)
267 (100%)
CM
55 (61%)
35 (39%)
90 (100%)
EM&CM
200 (63%)
119 (37%)
319 (100%)
NT
10 (71%)
4 (29%)
14 (100%)
NG
19 (49%)
20 (51%)
39 (100%)
NT&NG
4 (44%)
5 (56%)
9 (100%)
OS
0 (0%)
1 (100%)
1 (100%)
CERT
3 (100%)
0 (0%)
3 (100%)
AHO
1 (50%)
1 (50%)
2 (100%)
BI
0 (0%)
1 (100%)
1 (100%)
EA
79 (46%)
91 (54%)
170 (100%)
HO
12 (36%)
21 (64%)
33 (100%)
HB
4 (44%)
5 (56%)
9 (100%)
ICT
2 (29%)
5 (71%)
7 (100%)
MedV
3 (50%)
3 (50%)
6 (100%)
MobV
1 (50%)
1 (50%)
2 (100%)
TP
0 (0%)
3 (100%)
3 (100%)
TR
7 (70%)
3 (30%)
10 (100%)
TSL
7 (70%)
3 (30%)
10 (100%)
VNL
5 (100%)
0 (0%)
5 (100%)
ZW
32 (53%)
28 (47%)
60 (100%)
Vooropleiding
<0,001***
MBO
153 (48%)
165 (52%)
318 (100%)
HAVO
399 (56%)
309 (44%)
708 (100%)
VWO
26 (76%)
8 (24%)
34 (100%)
BD
2 (20%)
8 (80%)
10 (100%)
HO
16 (62%)
10 (38%)
26 (100%)
CD
32 (80%)
8 (20%)
40 (100%)
1 Mean (SD); n (%)
2 *p<0.05; **p<0.01; ***p<0.001
Toon code
## Laad dlookrsuppressMessages(library(dlookr))## Toon een samenvatting van de data, gesorteerd op missende waardendiagnose(dfOpleiding_inschrijvingen) |>mutate(missing_percent =round(missing_percent, 2),unique_rate =round(missing_percent, 2)) |>arrange(desc(missing_percent)) |> knitr::kable(caption ="Kwaliteit van de data voor bewerkingen (gesorteerd op missende waarden)",col.names =c("Variabelen","Type","# Missende waarden","% Missende waarden","# Unieke waarden","Ratio unieke waarden"))
Kwaliteit van de data voor bewerkingen (gesorteerd op missende waarden)
Uit de eerste diagnose blijkt dat niet alle variabelen goed genoeg zijn voor het bouwen van een prognosemodel: er zijn missende waarden en niet alle veldtypes zijn geschikt. We passen de variabelen aan zodat we in het model er goed mee kunnen werken.
Prognosemodellen kunnen niet omgaan met missende waarden. Om bias te voorkomen verwijderen we geen rijen met missende waarden, maar vullen die op (imputatie). We bewerken de data zo dat alle missende waarden worden opgevuld: bij numerieke waarden met het gemiddelde en bij categorische variabelen met ‘Onbekend’.
We passen sommige variabelen aan, zodat ze in het model gebruikt kunnen worden: tekstvelden zetten we om naar factor (een categorische variabele); logische variabelen (Ja/Nee) zetten we om naar een numerieke variabele (1/0).
De uitkomstvariabele, Retentie, leiden we af van de variabele SUC_Uitval_aantal_jaar_LTA. Als de waarde daar 1 is, is de student na 1 jaar uitgevallen, 2 na 2 jaar, etc. Zolang de waarde daar 0 is, is de student niet uitgevallen.
Een fictief studentnummer (INS_Student_UUID_opleiding_vorm) gebruiken we, zodat we - als er afwijkende resultaten zijn - de dataset gericht kunnen onderzoeken indien nodig.
Toon code
## Bewerk de datadfOpleiding_inschrijvingen <- dfOpleiding_inschrijvingen |>## Imputeer alle numerieke variabelen met de meanmutate(across(where(is.numeric), ~ifelse(is.na(.x),mean(.x, na.rm = T), .x )) ) |>## Zet character variabelen om naar factormutate(across(where(is.character), as.factor)) |>## Zet logische variabelen om naar 0 of 1mutate(across(where(is.logical), as.integer)) |>## Vul in factoren missende waarden op met "Onbekend"mutate(across(where(is.factor), ~suppressWarnings(fct_explicit_na(.x, na_level ="Onbekend") ))) |>## Herschik de kolommen, zodat Retentie vooraan staatselect(Retentie, everything()) ## Bekijk de data## glimpse(dfOpleiding_inschrijvingen) ## Laad dlookrsuppressMessages(library(dlookr))## Maak een diagnose van de datadiagnose(dfOpleiding_inschrijvingen) |>mutate(missing_percent =round(missing_percent, 2),unique_rate =round(unique_rate, 2)) |> knitr::kable(caption ="Kwaliteit van de data na bewerkingen",col.names =c("Variabelen","Type","# Missende waarden","% Missende waarden","# Unieke waarden","Ratio unieke waarden"))
Kwaliteit van de data na bewerkingen
Variabelen
Type
# Missende waarden
% Missende waarden
# Unieke waarden
Ratio unieke waarden
Retentie
factor
0
0
2
0.00
Aanmelding
numeric
0
0
282
0.25
Aansluiting
factor
0
0
6
0.01
APCG
factor
0
0
3
0.00
Cijfer_CE_Engels
numeric
0
0
59
0.05
Cijfer_CE_Engels_missing
factor
0
0
2
0.00
Cijfer_CE_Natuurkunde
numeric
0
0
22
0.02
Cijfer_CE_Natuurkunde_missing
factor
0
0
2
0.00
Cijfer_CE_Nederlands
numeric
0
0
48
0.04
Cijfer_CE_Nederlands_missing
factor
0
0
2
0.00
Cijfer_CE_VO
numeric
0
0
32
0.03
Cijfer_CE_VO_missing
factor
0
0
2
0.00
Cijfer_CE_Wiskunde
numeric
0
0
55
0.05
Cijfer_CE_Wiskunde_missing
factor
0
0
2
0.00
Cijfer_SE_VO
numeric
0
0
31
0.03
Cijfer_SE_VO_missing
factor
0
0
2
0.00
Collegejaar
numeric
0
0
11
0.01
Dubbele_studie
factor
0
0
2
0.00
Geslacht
factor
0
0
2
0.00
ID
factor
0
0
1135
1.00
Leeftijd
integer
0
0
19
0.02
Reistijd
numeric
0
0
475
0.42
SES_Arbeid
numeric
0
0
298
0.26
SES_Totaal
numeric
0
0
551
0.49
SES_Welvaart
numeric
0
0
391
0.34
Studiekeuzeprofiel
factor
0
0
22
0.02
Vooropleiding
factor
0
0
6
0.01
Toon code
detach("package:dlookr", unload =TRUE)
1.2.4 Bekijk de onderlinge correlaties
Het is verstandig om voorafgaand aan het bouwen van een model te kijken naar de onderlinge correlaties tussen numerieke variabelen. Dit geeft inzicht in de data en kan helpen bij het maken van keuzes voor het model of de duiding van de uitkomsten.
Toon code
## Maak een plot van de onderlinge correlaties in numerieke variabelendfOpleiding_inschrijvingen |>select(-Collegejaar) |>select(where(is.numeric)) |>cor() |> corrplot::corrplot(order ='hclust', addrect =4,method ="number", tl.cex =0.8, tl.col ="black",diag =FALSE)
1.2.5 Bouw de trainingset, validatieset en testset
De data is nu geschikt om een prognosemodel mee te bouwen.
Om het model te bouwen, testen en valideren, splitsen we de data in drie delen van 60%, 20% en 20%. We doen dit op zo’n manier, dat elk deel ongeveer een gelijk aantal studenten bevat dat doorstudeert (dus niet uitvalt).
We trainen het model op basis van 60% en valideren de modellen tijdens het trainen op de overige 20% (de validatieset).
De verdeling van de training- en validatieset muteren we 10x (10 folds) om te voorkomen dat het model te veel leert van de trainingset en daardoor slecht presteert op de validatieset.
Als het model klaar is, testen we het op de 20% studenten uit de testset. De testset blijft dus de gehele tijd ongemoeid, zodat we overfitting - een te goed model op bekende data, maar slechte presetaties (performance) op onbekende data - voorkomen.
Een willekeurig, maar vaststaand seed-getal voorkomt dat we bij elke run van het model c.q. deze code een net iets andere uitkomst krijgen.
Toon code
set.seed(0821)## Splits de data in 3 delen: 60%, 20% en 20%splits <-initial_validation_split(dfOpleiding_inschrijvingen,strata = Retentie,prop =c(0.6, 0.2))## Maak drie sets: een trainingset, een testset en een validatiesetdfRetentie_train <-training(splits)dfRetentie_test <-testing(splits)dfRetentie_validation <-validation_set(splits)## Maak een resample set op basis van 10 folds (default)dfRetentie_resamples <-vfold_cv(dfRetentie_train, strata = Retentie)
Verhouding training- en testset
Naam
Retentie
Aantal
Proportie
Trainingset
FALSE
304
44.7%
Trainingset
TRUE
376
55.3%
Testset
FALSE
102
44.7%
Testset
TRUE
126
55.3%
1.3 Model I: Logistische Regressie
Het eerste model is een logistische regressie met penalized likelihood; we gebruiken de glmnet engine voor het bouwen van het model. Penalized likelihood is een techniek die helpt bij het voorkomen van overfitting. Glmnet is een populair package voor het bouwen van logistische regressiemodellen.
We gebruiken de Area under the ROC Curve (AUC/ROC) als performance metric.
Vervolgens zetten we meerdere stappen in een ‘recipe’:
We definiëren de student-ID als ID variabele. Daarmee krijgt deze variabele de rol van uniek rij-kenmerk.
We verwijderen vervolgens de oorspronkelijke student-ID en het collegejaar uit de data, omdat deze verder niet gebruikt moeten worden in het model.
We converteren factoren naar dummy variabelen.
We verwijderen variabelen die geen waarde toevoegen: variabelen met enkel nullen.
We transformeren numerieke variabelen om ze met elkaar te kunnen vergelijken door ze te centreren en schalen.
Sterk gecorreleerde waarden verwijderen we nu niet, omdat we later in de analyse de eventuele samenhang met andere variabelen in een prognosemodel nog willen kunnen visualiseren.
## Bouw de recipe: logistische regressielr_recipe <-recipe(Retentie ~ ., data = dfRetentie_train) |>update_role(ID, new_role ="ID") |>## Zet de student ID als ID variabelestep_rm(ID, Collegejaar) |>## Verwijder ID en collegejaar uit het modelstep_dummy(all_nominal_predictors()) |>## Maak dummy variabelen van categorische variabelenstep_zv(all_predictors()) |>## Verwijder zero valuesstep_normalize(all_numeric_predictors()) ## Centreer en schaal numerieke variabelen## Toon de recipetidy(lr_recipe) |> knitr::kable(col.names =c("Nummer", "Operatie", "Type","Getraind","Sla over","ID"))
Nummer
Operatie
Type
Getraind
Sla over
ID
1
step
rm
FALSE
FALSE
rm_RuiEN
2
step
dummy
FALSE
FALSE
dummy_Nog78
3
step
zv
FALSE
FALSE
zv_o5wII
4
step
normalize
FALSE
FALSE
normalize_u1KRn
De variabelen die nu nog resteren zijn:
Resterende variabelen na bewerkingen
Aanmelding
APCG_Nee
Studiekeuzeprofiel_HO
Cijfer_CE_Engels
APCG_Onbekend
Studiekeuzeprofiel_HB
Cijfer_CE_Natuurkunde
Cijfer_CE_Engels_missing_Nee
Studiekeuzeprofiel_ICT
Cijfer_CE_Nederlands
Cijfer_CE_Natuurkunde_missing_Nee
Studiekeuzeprofiel_MedV
Cijfer_CE_VO
Cijfer_CE_Nederlands_missing_Nee
Studiekeuzeprofiel_MobV
Cijfer_CE_Wiskunde
Cijfer_CE_VO_missing_Nee
Studiekeuzeprofiel_TP
Cijfer_SE_VO
Cijfer_CE_Wiskunde_missing_Nee
Studiekeuzeprofiel_TR
Leeftijd
Cijfer_SE_VO_missing_Nee
Studiekeuzeprofiel_TSL
Reistijd
Dubbele_studie_Nee
Studiekeuzeprofiel_VNL
SES_Arbeid
Geslacht_V
Studiekeuzeprofiel_ZW
SES_Totaal
Studiekeuzeprofiel_CM
Studiekeuzeprofiel_Onbekend
SES_Welvaart
Studiekeuzeprofiel_EM.CM
Vooropleiding_HAVO
Retentie
Studiekeuzeprofiel_NT
Vooropleiding_VWO
Aansluiting_Direct
Studiekeuzeprofiel_NG
Vooropleiding_BD
Aansluiting_Na.CD
Studiekeuzeprofiel_NT.NG
Vooropleiding_HO
Aansluiting_Switch.extern
Studiekeuzeprofiel_CERT
Vooropleiding_CD
Aansluiting_Switch.intern
Studiekeuzeprofiel_AHO
Aansluiting_Tussenjaar
Studiekeuzeprofiel_EA
1.3.3 Maak de workflow
Voor de uitvoering bouwen we een nieuwe workflow. Daaraan voegen we het model en de bewerkingen in de recipe toe.
## Maak de workflow: logistische regressielr_workflow <-workflow() |>## Maak een workflowadd_model(lr_mod) |>## Voeg het model toeadd_recipe(lr_recipe) ## Voeg de recipe toe## Toon de workflowlr_workflow
Het model moet getuned worden. Dit houdt in dat we de beste parameters voor het model moeten vinden. We maken een grid met verschillende penalty waarden. Daarmee kunnen we vervolgens het beste model selecteren met de hoogste ROC/AUC. We plotten de resultaten van de tuning, zodat we hieruit het beste model kunnen kiezen.
## Maak een grid: logistische regressielr_reg_grid <-tibble(penalty =10^seq(-4, -1, length.out =30))## Train en tune het model: logistische regressielr_res <- lr_workflow |>tune_grid(dfRetentie_validation,grid = lr_reg_grid,control =control_grid(save_pred =TRUE),metrics =metric_set(roc_auc))
Toon code
## Plot de resultaten + een rode verticale lijn voor de max AUClr_plot <- lr_res |>collect_metrics() |>ggplot(aes(x = penalty, y = mean)) +geom_point() +geom_line() +## Maak de schaal van de x-as logaritmischscale_x_log10(labels = scales::label_number()) +theme(axis.title.x =element_text(margin =margin(t =20)) ) +# Bepaal de titel, ondertitel en captionlabs(caption = sCaption,x ="Area under the ROC Curve",y ="Penalty" )## Voeg LTA elementen toe lr_plot <-Add_LTA_Theme_Elements(lr_plot, title_subtitle =FALSE)# Zoek de penalty waarde met de max AUCmax_auc_penalty <- lr_res |>collect_metrics() |>filter(mean ==max(mean)) |>pull(penalty)# Voeg de rode verticale lijn toe aan lr_plotlr_plot_plus <- lr_plot +geom_vline(xintercept = max_auc_penalty, color ="red")# Vind een mean voor de max AUC die hoger ismax_auc_mean <- lr_res |>collect_metrics() |>filter(mean ==max(mean)) |>pull(penalty)## Print de definitieve plotlr_plot_plus
1.3.5 Kies het beste model
We evalueren modellen met een zo hoog mogelijke Area under the ROC Curve (AUC/ROC) en een zo laag mogelijke penalty. Zo kunnen we uit de resultaten het beste model kiezen. Tot slot maken we een ROC curve om de prestaties van het model te visualiseren.
## Toon het beste modeltop_models <- lr_res |>show_best(metric ="roc_auc", n =10) |>mutate(mean =round(mean, 6)) |>arrange(penalty) top_models|> knitr::kable(col.names =c("Penalty", "Metriek", "Estimator","Gemiddelde","Aantal","SE","Configuratie"))
## Verzamel de predicties en evalueer het model (AUC/ROC): logistische regressielr_auc <- lr_res |>collect_predictions(parameters = lr_best) |>roc_curve(Retentie, .pred_FALSE) |>mutate(model ="Logistisch Regressie")## Plot de ROC curveGet_ROC_Plot(lr_auc, position =1)
## Bepaal de AUC van het beste modellr_auc_highest <- lr_res |>collect_predictions(parameters = lr_best) |>roc_auc(Retentie, .pred_FALSE)## Voeg de naam van het model en de AUC toe dfModel_resultsdfModel_results <- dfModel_results |>add_row(model ="Logistic Regression", auc = lr_auc_highest$.estimate)
1.4 Model II: Tree-based ensemble
Het tweede model is een random forest: een ensemble van beslisbomen (decision trees). Het is een krachtig model dat goed om kan gaan met complexe data en veel variabelen.
We gebruiken de ranger engine voor het bouwen van het model.
1.4.1 Bepaal het aantal PC-cores
Omdat een random forest model veel berekeningen vereist, willen we daarvoor alle computerkracht gebruiken die beschikbaar is. Het aantal CPU’s (cores) van de computer bepaalt hoe snel het model getraind kan worden. Deze informatie gebruiken we bij het bouwen van het model.
Toon code
## Bepaal het aantal corescores <- parallel::detectCores()
1.4.2 Maak het model
We bouwen eerst het model. We gebruiken de rand_forest functie om het model te bouwen. We tunen de mtry en min_n parameters. De mtry parameter bepaalt het aantal variabelen dat per boom wordt gebruikt. De min_n parameter bepaalt het minimum aantal observaties dat in een blad van de boom moet zitten. De functie tune() is hier nog een placeholder om de beste waarden voor deze parameters - die we later bepalen - daar in te stellen. We gebruiken 1.000 bomen c.q. versies van het model.
## Bouw het model: random forestrf_mod <-rand_forest(mtry =tune(), min_n =tune(), trees =1000) |>set_engine("ranger", num.threads = cores) |>set_mode("classification")
1.4.3 Maak de recipe
We maken een recipe voor het random forest model. We verwijderen de student ID en het collegejaar uit de data, omdat deze niet moet worden gebruikt in het model. Overige stappen zijn bij een random forest minder relevant in tegenstelling tot een regressiemodel.
## Maak de recipe: random forestrf_recipe <-recipe(Retentie ~ ., data = dfRetentie_train) |>step_rm(ID, Collegejaar) ## Verwijder ID en Collegejaar uit het model## Toon de recipetidy(rf_recipe) |> knitr::kable(col.names =c("Nummer", "Operatie", "Type","Getraind","Sla over","ID"))
Nummer
Operatie
Type
Getraind
Sla over
ID
1
step
rm
FALSE
FALSE
rm_Q83Tz
1.4.4 Maak de workflow
We voegen het model en de recipe toe aan de workflow voor dit model.
## Maak de workflow: random forestrf_workflow <-workflow() |>add_model(rf_mod) |>add_recipe(rf_recipe)## Toon de workflowrf_workflow
══ Workflow ════════════════════════════════════════════════════════════════════
Preprocessor: Recipe
Model: rand_forest()
── Preprocessor ────────────────────────────────────────────────────────────────
1 Recipe Step
• step_rm()
── Model ───────────────────────────────────────────────────────────────────────
Random Forest Model Specification (classification)
Main Arguments:
mtry = tune()
trees = 1000
min_n = tune()
Engine-Specific Arguments:
num.threads = cores
Computational engine: ranger
1.4.5 Tune en train het model
We trainen en tunen het model in de workflow. We maken een grid met verschillende waarden voor de parameters mtry en min_n. We gebruiken de Area under the ROC Curve (AUC/ROC) als performance metric. Met de resultaten van de tuning kiezen we het beste model.
## Toon de parameters die getuned kunnen wordenrf_mod
Random Forest Model Specification (classification)
Main Arguments:
mtry = tune()
trees = 1000
min_n = tune()
Engine-Specific Arguments:
num.threads = cores
Computational engine: ranger
## Extraheer de parameters die getuned wordenextract_parameter_set_dials(rf_mod)
Collection of 2 parameters for tuning
identifier type object
mtry mtry nparam[?]
min_n min_n nparam[+]
Model parameters needing finalization:
# Randomly Selected Predictors ('mtry')
See `?dials::finalize` or `?dials::update.parameters` for more information.
## Bepaal de seedset.seed(2904)## Bouw het grid: random forestrf_res <- rf_workflow |>tune_grid(dfRetentie_validation,grid =25,control =control_grid(save_pred =TRUE),metrics =metric_set(roc_auc))
i Creating pre-processing data to finalize unknown parameter: mtry
1.4.6 Kies het beste model
We evalueren de beste modellen en maken een ROC curve om de performance van het model te visualiseren. Vervolgens vergelijken we de prestaties van de modellen en kiezen we het beste model.
## Toon de beste modellenrf_res |>show_best(metric ="roc_auc", n =15) |>mutate(mean =round(mean, 6)) |> knitr::kable(col.names =c("Mtry", "Min. aantal", "Metriek","Estimator","Gemiddelde","Aantal","SE","Configuratie"))
Mtry
Min. aantal
Metriek
Estimator
Gemiddelde
Aantal
SE
Configuratie
2
32
roc_auc
binary
0.620682
1
NA
Preprocessor1_Model03
2
6
roc_auc
binary
0.612200
1
NA
Preprocessor1_Model22
3
8
roc_auc
binary
0.596716
1
NA
Preprocessor1_Model08
4
26
roc_auc
binary
0.590181
1
NA
Preprocessor1_Model23
7
35
roc_auc
binary
0.589636
1
NA
Preprocessor1_Model14
6
23
roc_auc
binary
0.588702
1
NA
Preprocessor1_Model17
5
36
roc_auc
binary
0.586601
1
NA
Preprocessor1_Model09
8
31
roc_auc
binary
0.578431
1
NA
Preprocessor1_Model24
11
29
roc_auc
binary
0.578120
1
NA
Preprocessor1_Model02
12
13
roc_auc
binary
0.576797
1
NA
Preprocessor1_Model19
11
21
roc_auc
binary
0.576564
1
NA
Preprocessor1_Model15
13
25
roc_auc
binary
0.572985
1
NA
Preprocessor1_Model25
9
12
roc_auc
binary
0.572829
1
NA
Preprocessor1_Model05
15
20
roc_auc
binary
0.570339
1
NA
Preprocessor1_Model20
20
10
roc_auc
binary
0.567383
1
NA
Preprocessor1_Model18
## Plot de resultatenautoplot <-autoplot(rf_res) +theme_minimal() +labs(y ="roc/auc",caption = sCaption )## Voeg LTA elementen toe autoplot <-Add_LTA_Theme_Elements(autoplot, title_subtitle =FALSE)print(autoplot)
## Selecteer het beste modelrf_best <- rf_res |>select_best(metric ="roc_auc")rf_best|> knitr::kable(col.names =c("Mtry", "Min. aantal", "Configuratie"))
## Bepaal de AUC/ROC curverf_auc <- rf_res |>collect_predictions(parameters = rf_best) |>roc_curve(Retentie, .pred_FALSE) |>mutate(model ="Random Forest")## Plot de ROC curveGet_ROC_Plot(rf_auc, position =2)
Toon code
## Bepaal de AUC van het beste modelrf_auc_highest <- rf_res |>collect_predictions(parameters = rf_best) |>roc_auc(Retentie, .pred_FALSE)## Voeg de naam van het model en de AUC toe dfModel_resultsdfModel_results <- dfModel_results |>add_row(model ="Random Forest", auc = rf_auc_highest$.estimate)
1.5 De uiteindelijke fit
In de laatste stap van deze analyse maken we het model definitief.
We testen het model op de testset en evalueren het model met metrieken en de Variable Importance Factor (VIF).
1.5.1 Combineer de AUC/ROC curves en kies het beste model
Eerst combineren we de AUC/ROC curves van de modellen om ze te vergelijken. We kiezen het beste model op basis van de hoogste AUC/ROC.
Toon code
## Combineer de AUC/ROC curves om de modellen te vergelijkenGet_ROC_Plot(list(lr_auc, rf_auc))
Toon code
## Bepaal welke van de modellen het beste is op basis van de hoogste AUC/ROCdfModel_results <- dfModel_results |>mutate(number =row_number()) |>mutate(best =ifelse(auc ==max(auc), TRUE, FALSE)) |>arrange(number)## Bepaal het beste modelsBest_model <- dfModel_results$model[dfModel_results$best ==TRUE]sBest_model_auc <-round(dfModel_results$auc[dfModel_results$best ==TRUE], 4)
Het beste model is het Random Forest model met een AUC/ROC van 0.6207. Het Logistic Regression model heeft een AUC van 0.6144. Het Random Forest model heeft een AUC van 0.6207. We ronden de analyse verder af met het Random Forest model.
1.5.2 Maak het finale model
We maken het finale model op basis van de beste parameters die we hebben gevonden. Door in de engine bij importance de impurity op te geven, wordt het beste random forest model gekozen om de data definitief mee te classificeren.
## Test het ontwikkelde model op de testset## Bepaal de optimale parameters## Bouw de laatste modellenlast_lr_mod <-logistic_reg(penalty = lr_best$penalty,mixture =1) |>set_engine("glmnet") |>set_mode("classification")last_rf_mod <-rand_forest(mtry = rf_best$mtry,min_n = rf_best$min_n,trees =1000) |>set_engine("ranger", num.threads = cores, importance ="impurity") |>set_mode("classification")
1.5.3 Maak de workflow
We voegen het model toe aan de workflow en updaten de workflow met het finale model.
We voeren de finale fit uit. De functie last_fit past het model toe op de validatieset.
## Voer de laatste fit uitset.seed(2904)## Maak voor beide modellen een laatste fit, zodat we deze kunnen opslaan voor later gebruiklast_fit_lr <- last_lr_workflow |>last_fit(splits)last_fit_rf <- last_rf_workflow |>last_fit(splits)lLast_fits <-list(last_fit_lr, last_fit_rf) |>set_names(c("Logistic Regression", "Random Forest"))## Bepaal welk model het beste isif(sBest_model =="Logistic Regression") { last_fit <- last_fit_lr} elseif(sBest_model =="Random Forest") { last_fit <- last_fit_rf}## Bewaar de resultaten, de modelresultaten en de bijbehorende datasFittedmodels_outputpath <-Get_Model_Outputpath(mode ="last-fits")saveRDS(lLast_fits, file = sFittedmodels_outputpath)sModelresults_outputpath <-Get_Model_Outputpath(mode ="modelresults")saveRDS(dfModel_results, file = sModelresults_outputpath)sData_outputpath <-Get_Model_Outputpath(mode ="data")saveRDS(dfOpleiding_inschrijvingen, file = sData_outputpath)
1.5.5 Evalueer het finale model: metrieken en vif
We evalueren het finale model op basis van 4 metrieken: 1) accuraatheid, 2) ROC/AUC en 3) de Brier score (de Mean Squared Error) en 4) de Variable Importance Factor (VIF). Uit de VIF is op te maken welke variabelen het meest bijdragen aan de voorspelling van de uitkomstvariabele.
# Extraheer de feature importancedfVif <- last_fit |>extract_fit_parsnip() |> vip::vi() |>arrange(desc(Importance)) |>head(20)# Maak de plot met fill op de variabele 'Importance'importance_plot <- dfVif |>ggplot(aes(x =reorder(Variable, Importance), y = Importance, fill = Importance)) +geom_col(show.legend =FALSE) +## Maak de titel en captionlabs(title ="Meest voorspellende factoren",subtitle ="Op basis van de Variable Importance Factor (VIF)",x =NULL,y ="VIF-score",caption = sCaption) +theme_minimal() +Set_LTA_Theme() +theme(axis.title.x =element_text(margin =margin(t =20)) ) +coord_flip()## Voeg LTA elementen toe importance_plot <-Add_LTA_Theme_Elements(importance_plot, title_subtitle =TRUE)# Toon de plotprint(importance_plot)
1.5.6 Plot de ROC curve
Tot slot maken we een ROC curve om de prestaties van het definitieve model te visualiseren. De Sensitivity (True Positive Rate) en Specificity (True Negative Rate) worden hierin uitgezet. De Area under the ROC Curve (AUC/ROC) geeft de prestaties van het model weer. Het model scoort beter naarmate de AUC/ROC dichter bij de 1 ligt, de linker bovenhoek. De linker bovenhoek houdt in dat alle prognoses exact overeenstemmen met de werkelijkheid. Een AUC/ROC van 0,5 betekent dat het model niet beter presteert dan een willekeurige voorspelling.
## Toon de roc curveauc_lf <- last_fit |>collect_predictions() |>roc_curve(Retentie, .pred_FALSE) |>mutate(model ="Last fit")Get_ROC_Plot(auc_lf, position =3)
1.6 Conclusies
1.6.1 Het beste prognosemodel voor deze opleiding
Het beste prognosemodel blijkt het Random Forest model te zijn.
Van de prognosemodellen die we hebben ontwikkeld om retentie na 1 jaar te voorspellen, had het Random Forest model de hoogste AUC/ROC waarde (0.6207).
1.6.2 Mate van accuraatheid en lift
Een prognosemodel moet minimaal beter presteren dan een base-model om waarde op accuraatheid toe te voegen. Het base-model neemt de grootste klasse van de gemiddelde retentie na 1 jaar van de afgelopen jaren als basis. Stel we zouden tegen alle studenten zeggen dat ze hun studie gaan halen, dan is de mate van accuratesse gelijk aan dit base-model. Dit base-model is dus altijd hoger dan de 50% lijn van de AUC/ROC curve, tenzij het base-model toevallig precies 50% is.
De mate van accuraatheid van de toepassing van het model is laag (58.77%).
Base-model: 55.28% – Voor deze opleiding bereken we het base-model als volgt. Van alle studenten studeerde 55.28% door. De grootste klasse (en de accuratesse) van het base-model is daarmee (100% - 55.28% = ) 55.28% die doorstudeerde.
Accuratesse prognose: 58.77% – Het model voorspelt Retentie na 1 jaar met een accuratesse van 58.77%.
Lift: 3.49% – Het model scoort in de huidige opbouw met een verschil van 3.49% (de lift) iets beter dan de accuraatheid van het base-model (55.28%).
De prestaties van het model kunnen we verder uitdrukken in een confusion matrix. Hierin zien we de voorspellingen van het model en de werkelijke uitkomsten. De matrix geeft inzicht in de mate van correcte en incorrecte voorspellingen. Ter illustratie werken we de matrix uit voor een voorspelling waarop een bindend studieadvies (BSA) gebaseerd zou kunnen zijn.
We passen de confusion matrix nu toe op het model dat als beste naar voren kwam. De accuraatheid van dit model is 58,8%. De accuraatheid van het model berekenen we door de som van de diagonaal te berekenen: het aandeel goed voorspelde uitkomsten, Retentie = Retentie (True Positive) en Geen retentie = Geen retentie (True Negative), af te zetten tegen het totaal aantal voorspellingen: 43,9% + 14,9% = 58,8%. (NB. De weergave in deze confusion matrix is diagonaal gespiegeld vergeleken met het voorbeeld.)
Naast de accuraatheid van het model is het ook belangrijk om te weten welke factoren het meest bijdragen aan de voorspelling van retentie na 1 jaar. Daarin gaat de vergelijking met de prestaties van het basemodel mank. Dat model geeft op geen enkele manier aan waarom een student een kans op succes heeft, anders dan - ‘dit is gebruikelijk in deze opleiding’.
Ongeacht de mate van accuraatheid, is het voor ons onderzoek naar kansengelijkheid essentieel om te weten welke factoren het meest bijdragen aan de voorspelling van retentie na 1 jaar. Het gaat erom dat we het belang van de factoren in de voorspellingen kunnen begrijpen en duiden. Machine Learning is hiervoor uitstekend geschikt, omdat het de mogelijkheid biedt om de belangrijkste factoren en hun invloed te leren kennen (Shmueli, 2010; Shmueli & Koppius, 2011).
1.7 Vervolgstappen: Factoranalyse
De volgende stap (stap 2) is een verdiepende analyse van de mate waarin de factoren die we gevonden hebben van invloed zijn op Retentie na 1 jaar. We kijken naar de rangorde, of ze die doorstudeerde verhogen of juist verlagen en hoe stabiel de factoren zijn als we in andere volgordes aan het model toevoegen. Om het concreet te maken zullen we het model toepassen op een aantal fictieve studenten, die we opbouwen uit de meeste voorkomende waarden in deze opleiding. Dit is het onderwerp van analyse 2: de Factoranalyse.
Literatuur
Barocas, S., Hardt, M., & Narayanan, A. (2023). Fairness and Machine Learning: Limitations and
Opportunities. fairmlbook.org. http://www.fairmlbook.org
Shmueli, G., & Koppius, O. (2011). Predictive
Analytics in Information Systems Research. MIS
Quarterly, 35(3), 553. https://doi.org/10.2307/23042796
Verantwoording
Deze analyse maakt deel uit van het onderzoek naar kansengelijkheid van het lectoraat Learning Technology & Analytics van De Haagse Hogeschool: No Fairness without Awareness | Het rapport is door het lectoraat ontwikkeld in Quarto 1.4.549. | Template versie: 0.9.1.9000